Phase noise of dispersion-managed solitons 
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We quantify noise-induced phase deviations of dispersion-managed solitons (DMS) in optical fiber commu- 
nications and femtosecond lasers. We first develop a perturbation theory for the dispersion-managed nonlin- 
ear Schrodinger equation (DMNLSE) in order to compute the noise-induced mean and variance of the soliton 
parameters. We then use the analytical results to guide importance-sampled Monte-Carlo simulations of the 
noise-driven DMNLSE. Comparison of these results with those from the original, un-averaged, governing equa- 
tions confirm the validity of the DMNLSE as a model for many dispersion-managed systems, and quantify the 
increased robustness of DMS with respect to noise-induced phase jitter. 
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The performance of many lightwave systems is ultimately 
limited by quantum noise. Scientifically and technologically 
important examples include optical fiber communication sys- 
tems and femtosecond (fs) lasers: the former are a key en- 
abling technology for the information age, while Ti:sapphire 
fs lasers have applications to optical atomic clocks. Estimat- 
ing the performance of these systems is a timely problem. Be- 
cause both kinds of systems are designed to operate with very 
high accuracies, however, failures result from the occurrence 
of unusually large deviations, which makes calculating error 
rates extremely difficult. Direct Monte-Carlo (MC) compu- 
tations of failure rates are impractical due to the exceeding 
number of samples necessary to obtain reliable estimates, and 
analytical predictions are impossible due to the scale and com- 
plexity of these systems. In particular, errors in both systems 
are often phase-sensitive, and both systems employ the tech- 
nique of dispersion management, whereby pulses propagate 
through a periodic concatenation of components with opposite 
signs of dispersion [1, 2]. The probability of rare events can 
often be efficiently calculated using importance sampling (IS), 
with which the noise is sampled from a biased distribution that 
makes the rare events occur more frequently. For IS to be suc- 
cessful, however, one must bias towards the most likely noise 
realizations that lead to the events of interest. For systems 
modeled by the nonlinear Schrodinger equation (NLSE), this 
is made possible using soliton perturbation theory (SPT) [3- 
6], but this tool is not available in dispersion-managed (DM) 
systems. Recently [1 1] we developed a perturbation theory for 
the dispersion-managed NLSE (DMNLSE) that governs the 
long-term dynamics of DM optical systems [8-10], and we 
performed ISMC simulations of the pulse amplitude and fre- 
quency. Here we employ this perturbation theory in order to 
compute noise-induced means and variances, and we develop 
IS for the DMNLSE by explicitly formulating and solving the 
optimal biasing problem. We then perform ISMC simulations 
of the pulse phase, where the choice of biasing is non-trivial. 
Finally, we compare these results to the original, un-averaged 
system as well as to systems with constant dispersion. 

Perturbations of dispersion-managed solitons. The prop- 
agation of optical pulses in dispersion-managed fiber commu- 



nication systems [1] and Ti: sapphire lasers [2] is described by 
an equation which we refer to as NLSE+DM: 

ij- + \ d(z/Za) ^+g{z/Za)\q\ 2 q = l'V(t,z) . (1) 

Here z is the propagation distance, t is the retarded time, 
q(t,z) is the slowly varying electric field envelope (rescaled 
to account for loss/amplification), d(z/z a ) is the local disper- 
sion, and g(z/z a ) describes the periodic power variation due 
to loss/amplification. The choice of d(z/z a ) is called a disper- 
sion map, and z a is the dispersion map period. The forcing is 
v(t,z) = L^Li v„(f) 8(z - nz a ) , where 8{z) is the Dirac delta 
and V„(f) is white Gaussian noise, satisfying E[v„(f)] = and 
E[v„(f)v*,(r')] = o5(t -t')5 nn i, where E[-] denotes ensem- 
ble average, the asterisk complex conjugation, 8„„i is the Kro- 
necker delta and a 1 is the noise variance. 

Once the compression/expansion of the pulse in each dis- 
persion map is properly factored out, the core pulse shape 
obeys the DMNLSE [8, 10]. Namely, to leading order 
we can approximate the solution of Eq. (1) as q(co,z,Q — 
e- iC ^ a2 / 2 u(co,z), where f(co) = fe- im f(t)6t is the Fourier 
transform of /(f) (all integrals are complete unless limits are 

given), and £ = z/z a . Here C(Q = z a /{ (<*(£') - d) dC 
where d is the average dispersion. The exponential factor 
in front of u((0,z) accounts for the rapid "breathing", while 
the slowly varying envelope u((0,z) satisfies the perturbed 
DMNLSE: 

du , -d 2 u 
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+ II u (t+t') u (t+t") u *t+t' +t")R(t' ,t") df'df" = iv(t,z) , (2) 

where the asterisk denotes complex conjugate, and for brevity 
Ut t ) = u(t,z), etc. The kernel R(t',t") quantifies the average 
nonlinearity over a dispersion map mitigated by dispersion 
management: R(t't") = [f e im ' t ' +im " t " r(a)'a)")da)'da)", where 

r(*)= f>e" c (^(C)dC.' 

The DMNLSE and its solutions depend on a parameter s, 
called the reduced map strength, which quantifies the size 
of the dispersion variations around their mean. Dispersion- 
managed solitons (DMS) are traveling- wave solutions of the 
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DMNLSE. If u (t,z;s) = e a2 ^ 2 f{t ;s) , then f{(o) satisfies a 
nonlinear integral equation which can be efficiently solved nu- 
merically [11]. The invariances of the DMNLSE then yield 
from u (t,z;s) a four-parameter family of DMS: 

«dm S (f,z;s) = e @{t < z) Af(A{t - T);A 2 s) , (3) 

where A and £2 are the DMS amplitude and frequency, 
®(t,z) = &(t - T) +4> is the local phase, and T and 4> are 
respectively the mean time and mean phase. In the unper- 
turbed case, the mean time and mean phase evolve according 
to T = £2 and 4> = (A 2 +£2 2 )/2. [Hereafter, the dot denotes 
differentiation with respect to z.] 

In the presence of perturbations, the DMS will evolve. 
If u(t,z) — «dms + w solves the noise-perturbed DMNLSE, 
w(t,z) satisfies the corresponding perturbed linearized 
DMNLSE. But part of the noise goes to change the soliton 
parameters. The noise-induced DMS parameter changes at 
each map period are found by demanding that w(t,z) remain 
small, and are written as Q(nz+) = Q{nz„ ) + AQ, where for 
Q = A,a, T it is AQ = {e i& y Q , V n {t)) / i^yq) , while [11] 

A<P = {e<%Mt))/{yvy*) + n{e i %,v n (t))/{y £1 ,yn)- 

(4) 

Here (f,g) = Re ff*(t)g(t)dt is the inner product, 
y A (t), . . . ,y®(t) are the neutral eigenmodes and generalized 
eigenmodes of the linearized DMNLSE, and y_ A {t), ■ ■ ■ ,y$(t) 
are the adjoint modes [11]: yj — —dU/d£ , y& = i^U, 

» = «!, yA = j(u + ^ + 2s d Jf), (5) 

while y^ = iy A , y T = -iyn/A, = -y T /A, andv^ = U. Here 

Z=t- T(z), and U(t,z) = u{t,z)e~ m is the DMS envelope. 
All of these results reduce to those arising from soliton per- 
turbation theory for the NLSE [4] when 5 = 0. 

Noise-induced parameter variances. When the perturba- 
tion in (2) represents noise, the above results yield a system 
of nonlinear stochastic differential equations (SDEs) for the 
evolution of the DMS parameters under the effect of noise: 

Q = v Q {z), *=i(A 2 + £2 2 )+£2v r (z) + v 4 ,(z), (6) 

for Q = A,£l,T, where the source terms are Vg(z) = 
(e ,e y ,5) /(yQ,yo) for all Q. We employ a continuum 
approximation of Eqs. (6), considering v(t,z) to be a 
zero-mean Gaussian white-noise process with autocorrela- 
tion E[S(t,z)S*(t,z)] = a 2 8(t - t')8(z - z'). The sources 
Va(z), • • • , V<j)(z) are then independent zero-mean white-noise 
processes, with autocorrelation E[Sq(z)Sqi(z')] — Oq8(z — 
z'), where Oq — G 2 \\y q\\ 2 / {y Q,yo) 2 ■ All of these variances 
depend on the soliton amplitude A as well as on the map 
strength s, and therefore on the propagation distance z. As 
a result, it is not possible to integrate Eqs. (6) in closed form, 
even in the case of constant dispersion. If the amplitude devi- 
ations are not large, one can approximate a\, . . . , cj, as con- 
stant. In this limit, Eqs. (6) can be integrated exactly, to give 



G(z) - Qo + W e (z) for Q=A,a while 

T(z) = T + f£2(z')dz' + W T (z) , (7a) 
o 

4>(z) = \ f (A 2 (z')+Q. 2 (z'))dz' + fn(z')S T (z)dz' + Wt,(z), 
(7b) 

where Wq(z) = J z SQ(z')dz is a zero-mean Wiener process 
with autocorrelation E[Wg(z)Wg/(z')] = C7g5gg/ min(z,z') . 
Unlike other soliton parameters, the mean value of the soli- 
ton phase is affected by the noise: 

E[*(L)] = i (A 2 + £2 2 ) L + I (a 2 + a 2 ) L 2 . (8) 

Stochastic calculus also yields the variances of the noise- 
perturbed output soliton parameters as var[A(L)] = djh, 
var[£2(L)] = cr^L, var[r(L)] = CJ^L + |cj^L 3 and 

var[*(L)] = (cr| + £2 2 (7^)L + a„cj|(7^L 2 

+ i(A 2 cT 2 + ^cT 2 )L 3 + l(a A 4 + 4)L 4 . (9) 

The cubic dependence on distance of the phase jitter due to 
the Kerr effect is the Gordon-Mollenauer jitter [12], but note 
that additional contributions are present. Remarkably, these 
results are formally identical to those for the NLSE [13]. The 
dependence of the variance on the soliton amplitude, how- 
ever, is dramatically different, due to the different dependence 
on A of the norms and inner products [11]. More importantly, 
these results are not enough to accurately estimate the occur- 
rence of those rare events in which the noise produces large 
phase deviations, because: (i) the prediction for the mean 
phase is inaccurate, as we show below; (ii) the knowledge of 
noise-induced means and variances is not enough to estimate 
behavior in the tails, because not all soliton parameters are 
Gaussian-distributed; (iii) even if the output probability den- 
sity functions (PDFs) were Gaussian, extrapolating the results 
to reach the distribution tails would magnify all uncertainties 
exponentially, thereby making any prediction meaningless. 

Most-likely noice-induced phase deviations. Even though 
perturbation theory is not enough by itself to predict failure 
rates, it provides a key tool to implement IS. To successfully 
apply IS, one must first find the most likely noise realiza- 
tion subject to the constraint of achieving a given parameter 
change. For additive white Gaussian noise, this problem is 
solved by minimizing the negative of the argument of the ex- 
ponential in the noise PDF, namely the integral / |v„(x)| 2 <fx, 
subject to the constraint AQ n = Agtarget- The solution is [11]: 

WO =AQ^ l&[ %{t)(y_ Q ,y Q )/\\yJ [ 2 . (10) 

To induce a larger than normal parameter change, one can 
then bias the noise by concentrating the MC samples around 
V„, pt(*). That is, V^biasedO) = V„ )0p tW + v n(t), where 
V« jO pt(0 is given above and v„(f ) is unbiased. 

Once the most likely noise realization that produces a given 
parameter change AQ n at each map period is known, one must 
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also find the most likely way to distribute a total parame- 
ter change Agtot at the output among all map periods. In 
principle, when seeking large phase changes, one must bias 
an appropriate combination of all linear modes. Among the 
terms in the right-hand-side of Eq. (6), however, changes in 
Q. 2 and (e~ l& y T ,S)Q. are second-order in the noise, while 
changes in A 2 are first-order in the noise, because £1 = 
while A ^ 0. We thus introduce the auxiliary quantity (z) 
such that d<j)/dz=A 2 /2 + v$(z) and 0(0) =<E> , and consider 
the optimal biasing problem for 0(z). In the continuum ap- 
proximation, the biasing function is then 

Ht,z)=Ay^(y A ,yJ/\\yJ\ 2 + (t-A 2 /2)y^,yJ/\\yJ 2 . 

(11) 

[The direct phase biasing is not given by 0z a , but rather by 
(0 — A 2 /2)z a -] Minimizing the sum of the norm of this bi- 
asing function over all amplifiers is equivalent to finding func- 
tions A (z) and 0(z) that minimize the functional 
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(0 -A 2 12) 



dz. (12) 



The Euler-Lagrange equations associated with J[A, 0] yield 

,2/0 - 2 (13a) 

(13b) 



2A-^+A : 



dA 



<j>-A 2 /2 = co 2 
1 



2 d K] 



dA 



-2cA = 0, 



where c is a Lagrange multiplier. The solution of the system 
composed of Eqs. (13), together with the boundary conditions 
A(0) =A , 0(0) =A(L) =0 and 0(L) = 0t a rget, determines 
the optimal amplitude and phase paths around which one must 
bias the ISMC simulations. (The condition A (L) = applies 
because amplitude changes at z = L do not produce phase 
changes.) This system can be integrated numerically using 
relaxation methods or numerical continuation software. Dif- 
ferent output phases can be targeted by solving the system for 
different values of c, which determines the amount of bias- 
ing being applied (c = yields no bias). Equations (11-13), 
which are new, reduce to known results in the case of constant 
dispersion [4, 14]. But unlike the constant-dispersion case 
(and unlike the case of time biasing), here the direct phase bi- 
asing is not constant in z. Physically, this is a consequence of 
the different way in which noise is translated into phase jitter 
in the DMNLSE by way of the linear modes. 

ISMC simulations. We now discuss importance-sampled 
MC (ISMC) simulations aimed at computing the PDF of the 
soliton phase at the output. To quantify larger-than-normal 
phase deviations, we perform the following steps at each map 
period: (i) recover the underlying DMS from the noisy sig- 
nal; (ii) obtain the linear modes and adjoint modes of the 
linearized DMNLSE around the given DMS; (iii) generate 
an unbiased noise realization, shift its mean with the appro- 
priately scaled adjoint modes and update the likelihood ra- 
tios [11]. We then add the noise to the pulse, propagate the 
noisy signal to the next map period, and repeat this process 
until the signal reaches the output. For each noise realization, 



the full DMNLSE is used to propagate the signal. (The lin- 
earized DMNLSE is only used to guide IS via its modes.) 
Even though the noise-induced DMS parameter changes at 
each map period are small, the accumulation of these changes 
often results in a significantly distorted output signal. 

We choose system parameters based on realistic values for 
optical fiber communications. Typical values of system pa- 
rameters for fs lasers can be obtained from Ref. [15]. We con- 
sider a piecewise constant dispersion map, with equal-length 
normal and anomalous dispersion sections and local disper- 
sion coefficients of 23.27 and — 22.97 ps 2 /km, resulting in an 
average dispersion of 0.15ps 2 /km. We set the unit time to 
17 ps, corresponding to d = 1 and s = 4, and we use the result- 
ing dispersion length of 1,923 km to normalize distances. We 
consider a transmission distance of 6000 km (or L = 3.1201) 
and amplifiers spaced 100 km apart, for a total of N a = 60 with 
dimensionless spacing z a = 0.052, and we set the map period 
to be aligned with them. We take a nonlinear coefficient of 
1.7 (W-km)~', a peak power of 3.51 mW, a loss coefficient of 
0.25dB/km, a spontaneous emission factor of 1.65, resulting 
in a dimensionless noise variance a 2 = 1.873 • 10~ 3 and an 
OSNR of 9.3 dB. We normalize pulse powers with the power 
needed to have g = 1, namely 3.51 mW. Finally, we take input 
pulses to have unit peak amplitude, resulting in A = 1 . 

Numerical results and discussion. Figure 1 shows the nu- 
merically reconstructed phase mean and variance as a func- 
tion of distance for the DMNLSE, as well as the correspond- 
ing values for constant-dispersion NLSE solitons with same 
mean, amplitude or energy as the DMS, plus the predictions 
of perturbation theory. Note that the DMS has the lowest vari- 
ance of all. (The constant-dispersion NLSE soliton with same 
width as the DMS has a much lower energy, which makes it 
much more susceptible to Gordon-Haus jitter.) Note that the 
means and variances of the numerically reconstructed phase 
depend dramatically on the particular definition of phase used 
in the simulations. Hence, consistency is crucial to ensure 
agreement between theory and simulations. Here, the phase 
of a noisy pulse is defined (both in the theory and in the nu- 
merics) as that of the underlying DMS (obtained as in [11]). 

A significant discrepancy is evident between analytical and 
numerical results for the mean phase. No satisfactory expla- 
nation currently exists for this effect, which also occurs for the 
NLSE [14]. It is likely to depend on a failure of SPT and/or 
from second-order effects. (Numerical results show that the 
discrepancy also depends on the computational noise band- 
width.) On the other hand, the analytical prediction for the 
variance agrees very well with the numerical results, both for 
NLSE and DMLNSE. 

Figure 2 shows the PDF of the DMS output phase as com- 
puted from ISMC simulations of the DMNLSE (2), standard 
MC simulations of the NLSE+DM (1), plus a Gaussian distri- 
bution with variance given by perturbation theory and a PDF 
obtained from direct ISMC simulations of the SDEs (6). The 
ISMC results collect samples generated with a few biasing tar- 
gets, using multiple IS [7] to properly combine the data. The 
PDFs from both the DMNLSE and the NLSE clearly deviate 
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FIG. 1 : (Color online) Mean (top) and variance (bottom) of the soli- 
ton phase as reconstructed with standard MC simulations. Thick red 
lines: DMNLSE; blue lines: constant-dispersion NLSE soliton with 
same amplitude as the DMS; magenta lines: NLSE soliton with same 
energy; orange lines: NLSE soliton with same width. The dashed 
black lines show the predictions from perturbation theory. 




FIG. 2: (Color online) PDF of output phase. Thick red curve: 
ISMC simulation of DMNLSE with 50,000 samples. Blue squares: 
standard MC simulation of the NLSE+DM with 250,000 samples. 
Dashed curve: ISMC simulation of the noise-driven SDEs (6). Dot- 
dashed curve: a Gaussian PDF with variance given by (9). 

from Gaussian, but they agree very well with each other as 
far down in probability as the unbiased MC simulations can 
reach. Conversely, while the Gaussian approximation agrees 
well near the peak of the PDF, for deviations from the mean 
phase of K or more (a value that is relevant for fs lasers) it is off 
by several orders of magnitude. Similarly, the SDEs obtained 
from perturbation theory fail to accurately reproduce the full 
dynamics of the soliton phase at lower-than-average values of 
phase. Remarkably, however, ISMC simulations guided by 
perturbation theory yield the correct phase behavior. 

Importantly, results from the noise-perturbed DMNLSE 



and NLSE+DM agree pathwise, not just in the overall PDFs 
at the output. That is, they agree for each noise realization 
as a function of distance. These results, which are surprising 
given the "softness" of the phase and the complexity of the 
system (nonlinearity, dispersion, noise, large deviations etc.), 
provide further confirmation of the validity and robustness of 
the DMNLSE in capturing the essential dynamics of DM sys- 
tems. Its usefulness is also increased by the availability of 
tools such as the perturbation theory presented here, an ana- 
logue of which is lacking for the NLSE+DM. 

To our knowledge, this is the first time that the probabil- 
ity of large phase deviations in dispersion-managed systems 
has been quantified. Similar dynamics should also arise for 
non-solitonic pulses, but the analysis for that case will be 
more complicated because generic pulses do not preserve a 
flat phase across their temporal profile upon propagation. 

An important question is also whether these results can be 
used in fs lasers in order to quantify the probability of the 
occurrence of phase slips in optical atomic clocks. Since gain 
and loss play an obvious role in lasers, one could expect that it 
will be necessary to derive a perturbation theory for the non- 
conservative version of the DMNLSE that was derived as a 
model for fs lasers [16]. Since the DMNLSE itself provides a 
surprisingly good quantitative description of these lasers [15], 
however, whether or not such an extension will indeed be nec- 
essary remains at present an open question. 
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